Critical Casimir forces and colloidal aggregation: A numerical study 
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We present a numerical study of the effective potential V e s between two hard-sphere colloids 
dispersed in a solvent of interacting particles, for several values of temperature and solvent density, 
approaching the solvent gas- liquid critical point. We investigate the stability of a system of particles 
interacting via V e s to evaluate the locus of colloidal aggregation in the solvent phase-diagram, and 
its dependence on the colloid size. We assess how the excluded volume depletion forces are modified 
by solvent attraction and discuss under which conditions solvent critical fluctuations, in the form of 
critical Casimir forces, can be used to effectively manipulate colloidal aggregation. 
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INTRODUCTION 

Understanding effective interactions between colloidal 
particles in the presence of macromolecular additives (e.g. 
surfactants, polymers, micelles, etc.) is of fundamental 
importance in soft-matter science. Since the pioneer- 
ing works of Asakura and Oosawa (AO) pQ and Vrij [2], 
who modeled the additives as an ideal gas depletant, it is 
known that the presence of (small) co-solutes in solution 
gives rise to entropic forces that significantly modify the 
structure and the dynamics of colloidal systems [3j [4] . 
These depletion forces arise from the entropy-driven ex- 
clusion of the co-solute from the volume between particles 
when their surface-to-surface relative distance is compa- 
rable to the co-solute size [5] . Accurate numerical [6j [7] 
and theoretical [8HTD] calculations of the effective poten- 
tial between hard-sphere colloids dispersed in a depletant 
of small hard-sphere (HS) particles have shown that the 
entropic attraction arising at short distances is followed 
by an oscillatory behavior due to excluded volume corre- 
lation effects. 

It is also well known that depletion forces can be strong 
enough to induce colloid aggregation. Decreasing the size 
ratio q between the hard-sphere depletant (of diameter 
a s ) and the hard-sphere colloid (a c ) gives rise to a critical 
point [7 at which colloids separate in a colloid-rich and 
a colloid-poor phase, the analog of the gas-liquid separa- 
tion in one-component systems. Interestingly, the smaller 
q = <7 s /<t c , the smaller is the volume fraction of small par- 
ticles (j) s = (tt /6)(N s /V)(Ts (where N s is the number of 
depletant particles and V is the available volume they can 
access) required to induce phase separation, bringing to 
the paradoxical predictions that it would seem impossi- 
ble to dissolve colloids in the presence of extremely small 
additives. Clearly, although modeling all interactions in 
terms of excluded volume provides a model amenable to 
theoretical and numerical treatment, this approach does 
not fully capture the ingredients behind effective poten- 
tials in real systems. This motivates recent attempts to 



include colloid-depletant [TT] and depletant-depletant in- 
teractions [T21IT5] . It has been shown that the latter sig- 
nificantly affect the resulting effective potential, altering 
the strength of the attraction and, in some cases, sup- 
pressing the oscillatory behavior characteristic of the de- 
pletant granularity [T2HT5] . 

A more realistic description of depletant-depletant in- 
teractions brings in two important features: first, beside 
S , the temperature T becomes a relevant variable for 
modulating the effective potential; second, the depletant 
itself can undergo a gas-liquid phase separation. Close 
to the critical point of the depletant, critical fluctua- 
tions —with the associated critical Casimir forces [T6l - 
[T8]— might become the dominant contribution in the ef- 
fective potential -with respect to the depletion force - 
hence driving colloidal aggregation [T9H2T] . For this to 
take place, as recently pointed out [22], it is necessary 
that the colloidal solution remains stable outside the de- 
pletant critical region, i.e. that Casimir-driven phase- 
separation is not preempted by other effects. A recent 
experimental work [23] has addressed the role of criti- 
cal fluctuations in colloidal aggregation. In this study 
colloidal particles were dissolved in a micellar solution 
which is known to undergo a gas-liquid phase separation 
upon heating. Exploring the onset of colloidal aggrega- 
tion as a function of T and micellar concentration, the 
authors were able to show that a continuous line can be 
drawn in the (0 S ,T) plane separating the region where 
colloids are stable from the region where colloids aggre- 
gate. This line was found to monotonically connect the 
pure depletion limit to the micellar critical point, sug- 
gesting the possibility that critical fluctuations can be 
considered as an extreme case of depletion forces where 
the dense critical regions act as renormalized depletants, 
thus proposing a link between depletion forces and crit- 
ical Casimir effect [23]. Although this interpretation in 
terms of a unique description of forces originated by dif- 
ferent mechanisms is fascinating, critical Casimir forces 
have a richer phenomenology when compared to standard 
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depletion forces. Indeed critical Casimir forces depend 
strongly on the interaction between the colloid and the 
critical depletant and it is also possible to generate re- 
pulsive Casimir forces [18 j, a case which does not have a 
depletion analog. 

Numerical simulations of simple models can help un- 
derstanding if, how and when critical fluctuations may 
play a role in colloidal aggregation. To this aim, in this 
work we evaluate numerically the effective potential be- 
tween two hard-sphere colloids immersed in an implicit 
solvent, in the presence of interacting depletant parti- 
cles, for a wide region of state points (0 S ,T), including 
the critical region. We investigate two different depletant 
models at several values of q, between 0.05 < q < 0.2, at 
the limit of todays computational possibilities, in order 
to disentangle aggregation driven by critical fluctuations 
from that arising from different phenomena. We provide 
a complete evaluation of the colloid stable and unstable 
regions and the relative location of these regions with re- 
spect to the depletant gas-liquid critical point as a func- 
tion of q. Our main result is that critical fluctuations 
are responsible for colloidal aggregation only under very 
specific conditions. 



MODELS AND METHODS 

We calculate the effective forces between two hard- 
spheres immersed in two different models for the deple- 
tant particles. In the first model, depletant particles in- 
teract via pairwise square- well potential (SW) 



oo, r < a s 

-e, a s <r < (1 + A)cr s (1) 
0, r>(l + A)a s 



Tij/\rij\^ where is the vector connecting the centers 
of two particles. The value of cos# = 0.8947 controls 
the width of the patch, while the value of A = 0.119 
fixes the interaction range. Therefore two patchy parti- 
cles interact via an attractive SW potential only if the 
two patches are properly aligned. The depletant parti- 
cle diameter a s and the potential well e (equal for both 
models) are chosen as unit of length and energy. T is 
measured in units of e (i.e. Boltzmann constant ks — 1). 
Both models reduce to the HS model at high T, when 
attraction becomes irrelevant. The SW and the 3P mod- 
els are characterized by a gas-liquid critical point, lo- 
cated respectively at (T c = 0.478, <j> c = 0.25 [26 ) and 
(T c = O.1247,0 C = 0.073). The difference in the critical 
<j) between the two models, due to the very small (j) c of 
the 3P model, provides a different behavior for the sol- 
ubility of colloids along the critical isochore of the two 
depletants. The investigation of both models thus offers 
the possibility to establish this behavior for the case of 
a simple depletant (SW) as well as for a depletant char- 
acterized by a small value of (j) c as in the experimental 
system of Ref. 23. We perform Monte Carlo simulations 
of two colloids fixed at several equally-spaced distances 
in a depletant of SW or 3P particles, as a function of 
T and depletant (reservoir) concentration. The size of 
the rectangular box has been chosen in such a way that 
the depletant density far from the colloids reaches the 
constant reservoir (bulk) value, in all directions. The 
effective force F e ^ is evaluated by performing virtual dis- 
placements of Ar (both positive and negative) for each 
colloid and computing the probability of having at least 
one collision P over i ap (Ar) = (1 — e - ^ hs ) with the deple- 
tant particles. As shown in Ref. [27 : 



with a well width A = 0.1 and a well depth e. The sec- 
ond depletant model is an anisotropic three-patches (3P) 
Kern-Frenkel 25] system which consist of hard-sphere 
particles decorated with three attractive sites whose in- 
teraction potential can be defined as: 

(rij-hi.fij) = ^(nj) • /(%n,,n,) (2) 

i.e. as a SW potential modulated by an angular function 

fij • n- a) > cos(#) 



1, if- 



and —Tij • fr- > cos(#) 
otherwise. 



(3) 

In EqJ3J is the unit vector associated to the orien- 
tation of patch a = 1,2,3 on particle and f^ = 



/3Feff(r) 



lim 

Ar^0+ 



r i ap (Ar) 



Ar 



lim 



r i ap (Ar) 



Ar 



(4) 



We have chosen Ar in such a way that P over i ap (Ar) ~ 5%. 
The effective colloid-colloid potential /3V e s (r) is then ob- 
tained by integrating /3F G s(r) from oo to r. In addi- 
tion, we perform standard MC simulations of the bulk 
SW model to numerically study the critical region. Also, 
to assess colloidal stability upon varying depletant con- 
ditions, we perform grand canonical Monte Carlo simula- 
tions of colloidal particles interacting via V e f / for several 
values of q and <\> c . 
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RESULTS 



compressibility) 



Critical behavior of the depletant 

Before discussing the results for the effective forces, we 
investigate the critical region of the pure depletant solu- 
tion. To quantify the growth of the thermal correlation 
length £ on approaching the critical point we perform sim- 
ulations of a bulk SW model along the critical isochore [26 
for several T values and evaluate numerically the spheri- 
cally averaged structure factor S(k) [35], shown in Fig. [I] 
The correlation length £ is estimated fitting the low wave 
vector k behavior of S(k) according to 



S(k) 



S(0) 



1 + 



(5) 



In Eq. [5j a constant c is added to the critical Lorentzian 
shape to account for the low q non-critical component, 
which is still not negligible at the highest investigated T. 
The best fitting functions are reported also in Fig. [I] 
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FIG. 1. Static structure factors S(k) (symbols) for a bulk 
SW system of 40000 particles along the critical isochore as a 
function of T. Dashed lines are fit according to Eq|5] The 
fitting parameter c of Eq[5] is found to be constant at all T 
and equal to c c± 0.22. The fit is restricted to a fixed interval 
(for all T) delimited by the two vertical dashed lines (0.4 < 
ka s < 3.5). The small k values of S(k) are indeed affected by 
numerical errors related to the small number of independent 
wave vectors contributing to the spherical average. The large 
k cut-off excludes the non-critical part of S(k). 



In the critical region a universal power-law behavior is 
expected for £ and 5(0) = pksTxr (where xt is the 



S(0) * Co 



T-T c 
T-T r 



(6) 
(7) 



with v = 0.63 and 7 = 1.24, since the SW model belongs 
to the Ising universality class. £0 and Co are the non- 
universal critical amplitudes. 
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FIG. 2. Power- law behavior of the thermal correlation length 
(top panel) and S(0) (lower panel) approaching the critical 
point for a bulk SW system (circles), resulting from fitting 
S(k) according to Eq. [5] Lines are power-law fits, with best 
fit exponents v — 0.65 and 7 = 1.20, in good agreement with 
the expected Ising values. The best fit value for the ampli- 
tude £0 ~ 0.23cr s . Squares indicates the characteristic length 
scale obtained from fitting the effective potential according to 
Eq.[l0] 

The two panels in Fig. [2] show that both £ and 5(0) 
indeed satisfy Eq. [6] with v and 7 values consistent with 
the expected Ising ones, when (T — T c )/T c < 0.1. Thus, 
the onset of the critical region can be defined as T w 
\.\T r . 



Effective Potentials 

Fig. [3] shows the evolution of the effective potential 
along the reservoir critical isochore. As previously ob- 
served for the Baxter model [13], V e & progressively loses 
its oscillatory character upon cooling and gradually turns 
into a completely attractive potential. In addition, a 
progressive significant increase of the interaction range 
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FIG. 3. Effective potential for g = 0.1 for the (a) SW and 
(b) 3P potentials along their respective critical isochores for 
different values of T* = T/T c . Dashed lines for T* < 1.045 are 
fits according to Eq. [To] with fit parameters Asw = —0.63, 
A 3 p -0.19 and £sw 2.23, f 3 p 2.28. 



is observed on approaching T c , signaling the onset of 
critical Casimir forces. A continuous path thus connects 
the high-T pure HS depletion to the Casimir interaction. 
Casimir forces are very dependent on the geometry in 
which the critical fluid is confined as well on the bulk 
universality class of the (confined) fluid and the bound- 
ary conditions of the confining surfaces [T7J [36] . Recently 
it was shown that in the case of a critical binary mixture 
confined between two spherical surfaces of diameter a at 
minimum distance between the two surfaces z = r — a, 
the potential generated by the critical Casimir forces can 
be written as 



0*(z) = -€>(z/Z). 



(8) 



Eqj8] is valid only if a ^> z , i.e. if one can exploit the 
Derjaguin approximation [29 . 0(z/£) is a scaling func- 



tion which contains all the information about the critical 
properties of the system, the geometry and the boundary 
conditions, while £ is the thermal correlation length. The 
sign of the Critical Casimir forces (repulsive or attractive) 
depend on the so called boundary conditions, i.e. on the 
ability of the two large colloids of absorbing or not the 
depletant. In the present case, the interaction between 
the colloid and the SW or 3P particles is purely repulsive 
(hard-sphere). We indicate in the following such condi- 
tion with the symbol ( ). The long distance behavior 

of (sphere-sphere) 0(z/£)(__) (i.e. for z/£ ^> 1) is [17] 



e(z/0(—)(z/^l) = 7TA(z/0e- 



(2/0 



(9) 



where A is a universal constant. It follows that the long 
distance behavior of the critical Casimir potential be- 
tween the two spherical colloids behaves as 



(7c, 



/3$(r) = TrAy e 



<r-* c )/£ 



(10) 



Fig. [3] shows that close to criticality, the effective po- 
tential /3V e ff for both the SW and 3P models can be 
represented by EqjTol in the entire range of distances. 
Fitting the low T numerical results for fiV e ff with Eq 10 
it is possible to extract a value of £ that should corre- 
spond to the correlation length of the critical bulk fluid. 
To show that the critical Casimir forces are controlling 
the effective interactions near T c we compare in Fig. [2] 
the £ value resulting from the fit of the effective potentials 
with that of the bulk SW system. The good agreement 
between the two values of £ strongly supports the idea 
that the effective potential close to the critical point is 
controlled by critical fluctuations. 




FIG. 4. Effective potential for q = 0.1 for the 3P potential 
and the SW (inset) for three different values of the reservoir 
cj) Sl along an isotherm. 
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For completeness, Fig. [4]shows the dependence of V e ff on 
depletant concentration for the 3P model when q = 0.1. 
The absolute value of the effective potential at contact be- 
comes larger with increasing S , due to the larger amount 
of depletant in suspension. Instead, the interaction range 
has a maximum at the critical isochore, significantly in- 
creasing only on approaching T c , retaining the typical de- 
pletion interaction range when far from the critical point. 
The inset of Fig. [4] shows similar data for the SW model. 



Phase separation of colloids: role of the depletant 
size 

To evaluate the stability of the colloidal solution 
in different regions of the depletant (reservoir) phase- 
diagram we investigate a system of colloidal particles in- 
teracting via the pair-wise additive Hamiltonian H = 
Y.ij V Hs(rij) + V e ff(rij) (where nj is the relative dis- 
tance between two generic particles i and j). While 
the assumption of pairwise additivity is essentially un- 
controlled, since the presence of a third particle in the 
vicinity of a pair of colloids will alter the depletant spa- 
tial distribution, it is expected that such approximation 
becomes progressively less relevant on decreasing q and 
that its main effect is to slightly shift the coexistence 
line towards higher values of the depletant density in the 
reservoir [34] . We thus perform grand canonical MC sim- 
ulations (GCMC) for several values of the chemical po- 
tential, spanning the entire density region. In this way we 
are able to detect if particles interacting via /3Vhs + PVett 
form a stable fluid phase for all densities or if there exist a 
density region where particles phase separate in colloidal- 
rich colloidal-poor coexisting phases, the analog of the 
gas-liquid coexistence. Repeating such calculations for 
the effective potential evaluated at different (j) s and T 
state points, the depletant bulk phase diagram can be 
separated into a region in which colloids are in a homo- 
geneous stable fluid phase for all densities from a region 
in which colloids undergo phase separation. 

The loci separating fluid from phase separated state 
points, according to GCMC calculations, for different q 
values are displayed in Fig. [5] in the phase diagrams of 
both depletant models. In each figure the gas-liquid crit- 
ical point and the phase coexistence region are reported. 
Each locus starts from the HS pure depletion limit and 
converges toward the critical point, ending on the gas 
side of the depletant gas-liquid coexistence curve. In the 
case of the SW, for q = 0.1, the crossover T at which the 
solution starts to be unstable along the critical isochore 
is much higher than the T of the onset of critical regime. 
Colloidal phase separation in this case is therefore not 
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FIG. 5. Phase diagram of (a) the SW model and of (b) 
the 3P model in the T/T c and reservoir packing fraction (f) s 
plane. The diamond indicates the critical point. The area 
below the critical point (green) is the phase separation region. 
Curves interpolating circles, squares and triangles are the loci 
of colloidal aggregation for several q values. The q = 0.025 
case in (b) has been calculated assuming the validity of the 
Derjaguin approximation. The dashed lines connecting the 
points are guides to the eye. The reduced temperature at 
which the solvent second virial coefficient is zero is T* = 1.506 
for the SW model and T* = 1.738 for the 3P model. 



driven by depletant criticality but by standard depletion 
forces. In the case of the 3P depletant of Fig. [5] (b), the 
aggregation line for q = 0.1 intersect the gas-liquid coex- 
istence very close to the depletant critical point, suggest- 
ing that, at the depletant critical density, it is possible 
to go close enough in T to the critical point to induce 
colloid phase separation via Casimir forces. We remark 
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that such a difference originates from the different critical 
packing fraction between the two models. 

We also study how the locus of onset of instability de- 
pends on q. To this aim, we repeat our calculations for 
q = 0.05 (for SW also q = 0.2), and the corresponding re- 
sults are also shown in Fig. [5] At this small q value, more 
than 30.000 depletant particles are needed in the simula- 
tions, requiring a significant numerical investment. It is 
possible to extend these calculations to smaller q values 
(where simulations becomes prohibitive due to the large 
number of solutes) by invoking the Derjaguin approxi- 
mation [29] , according to which [37] for two different size 
ratios q\ and q2 holds: 



Ves(h/a a ,qi) V e ^(h/a s ,q 2 ) 



(?2 



(11) 
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IS 



where h is the surface to surface distance. Eq. 
expected to hold accurately as far as the curvature radius 
of the interacting spheres is sufficiently larger than the 
range of the interaction potential. 

We find (not shown) that when the ratio between the 
interaction range and the colloid diameter is smaller than 
10% the approximation is rather good. In principle, we 
can thus extend the q = 0.05 data to smaller q and draw 
stability lines for arbitrary values of q < 0.05, as shown 
for example in Fig. [5jb) for q = 0.025. Our results show 
that q plays a significant role in the location of the insta- 
bility line, in agreement with the results reported in a pi- 
oneering, but largely unnoticed, experimental study [30] . 

Finally we notice that along the stability lines the value 
of the reduced second virial coefficient, i.e. the second 
virial coefficient normalized to the hard-sphere value [28] , 



B* 



3_ r- 



[l-exp(-/3V eff )]r 2 dr 



(12) 



is close to the value B% 



-1.6 which has been identified 



as the universal value for the onset of gas-liquid separa- 
tion in short-range attractive potentials [28] . 



DISCUSSION AND CONCLUSIONS 

Several considerations can be drawn from the analy- 
sis of the data shown in the previous figures, (i) The 
depletant-depletant attraction always favors colloidal ag- 
gregation as shown by the monotonic behavior of the 
curves in Fig. [5] Along the critical isochore, approach- 
ing the critical T, /3V G s becomes increasingly long-ranged 
and its radial dependence is well represented by an ex- 
ponential decay. This decay length is controlled by £, in 



agreement with theoretical predictions (Eq. 10 ) [17J . The 



amplitude is found to be T-dependent stressing the im- 
portance of the condition £ ^> a s for fully entering in 
the theoretically expected scaling regime, (ii) /3V e s has 
been used to evaluate the colloidal aggregation loci under 
the pairwise additivity approximation. We find that the 
g-dependence of these lines is quite strong. Indeed, each 
locus starts from the high T HS limit, which provides a 
lower value of <p s for smaller q [7], and bends towards 
lower (j) s until merging with the depletant gas-liquid co- 
existence, in agreement with the experimental work of 
Ref. [23 j. The marked dependence of the locus of ag- 
gregation on the size of the colloid strongly suggests that 
no universal connection can be exploited between spe- 
cific thermodynamic properties of the depletant (which 
are <j c independent) and aggregation. Interestingly, for 
the SW depletant, colloids as small as ten times the de- 
pletant size (see q = 0.1 or smaller in Fig. [5^i), aggregate 
well before the depletant develops significant critical fluc- 
tuations, proving that colloidal aggregation driven by a 
critical depletant requires specific conditions to be ob- 
served [22]. However, depletants with a low <\> c (as exhib- 
ited by particles with limited valence interactions [31] [32]) 
enhance the possibility of observing Casimir aggregation 
since the underlaying depletion interaction, whose am- 
plitude is controlled by <fr s , is weaker. This is indeed 
the case for the 3P model for q = 0.1 (see Figj5]o) as 
well as for the experimental system investigated in Ref. 
23. (Hi) Finally, we observe that critical Casimir forces 
may play a dominant role in colloidal aggregation when a 
colloid-colloid repulsive interaction compensates the non- 
critical component of the depletant mediated attraction 
(which is limited to the <j s scale). Indeed, the effective 
potential calculated here can also be used to evaluate 
the stability line when the colloid-colloid hard-sphere in- 
teraction is complemented by additional terms. When 
an additional repulsive potential of the order of a few 
/c#T, acting on the length scale of a few depletant di- 
ameters, is present (the typical case for screened elec- 
trostatic interactions), we find that colloidal aggregation 
takes place only very close to the depletant critical point 
(i.e. only when £ ^> a s ), where the critical Casimir poten- 
tial, which decays on the scale of £, becomes the leading 
term in /3V e fi. Hence, our study goes beyond the purpose 
of using colloids as model systems for studying Casimir 
forces [18, 33 since it suggests that with an appropriate 
tuning of electrostatic interactions and depletant condi- 
tions, critical Casimir forces can be exploited to finely 
control collective assembly [2T| [23] . The sensitivity of the 
depletant critical fluctuations, where strength and range 
of interaction can be significantly changed via a tiny T 
variation, opens up a new way of manipulating colloidal 
aggregates which need to be exploited in the near future. 
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